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Abstract 

Although Potent purports to use only radial velocities in retrieving the potential ve- 
locity field of galaxies, the derivation of transverse components is implicit in the smoothing 
procedures. Thus the possibility of using nonradial line integrals to derive the velocity 
field arises. In the case of inhomogeneous distributions of galaxies, the optimal path for 
integration need not be radial, and can be obtained by using max-fiow algorithms. In this 
paper we present the results of using Dijkstra's algorithm to obtain this optimal path and 
velocity field. 

1 Introduction 

One of the most effective methods for reconstructing peculiar velocity fields of galaxies on the 
scale of 100 Mpc from observed redshifts has been POTENT ([1], [3] hereafter DBF, [4]). The 
main weakness of this method is the sparseness of the data. Sky coverage of the most complete 
redshift surveys is still highly anisotropic. Some of this incompleteness will be overcome by 
more rigorous and deeper surveys. However the intrinsic distribution, clustering of galaxies 
and voids will always be a barrier to determining peculiar velocity fields accurately. In these 
circumstances it is important to find the optimal way of reconstructing the potential peculiar 
velocity field. This in turn could lead us to a much improved understanding of local density 
inhomogeneities, and further constraints on the value of and on the various cold and hot 
dark matter model candidates. 

Potent has usually been presented as a method for obtaining the smoothed peculiar ve- 
locity field almost directly from observations of redshifts and distances of galaxies. Such an 
attractive idea often obscures the actual practice by which the smoothed potential field is con- 
structed, and in so doing obscures also ways in which the method can be improved. In this 
paper we present a method for obtaining the velocity potential from non-radial paths. This 
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Figure 1: The optimal path minimises the error of the velocity potential, and will be pulled 
towards regions of high galaxy number density. 

opens up the possibility of choosing paths that pass through regions where the number density 
of galaxies is high, and where the errors of distance estimation are low. Such paths should yield 
a more accurate potential velocity field. To understand this new approach, we shall have to 
outline the main steps of the POTENT method of DBF. 

2 Orthodox Potent 

The main steps of POTENT may be crudely summarised as the following 

1. Take measured redshifts and distances of galaxies and hence obtain peculiar velocities of 

individual galaxies. 

2. Construct an initial smoothed velocity field by using a tensor window function. 

3. Take the radial component of this initial smoothed velocity field and carry out the line 

integral along a radial path to obtain the velocity potential. 

4. Use the potential to derive the vectorial smoothed peculiar velocity field. 

Of course there are numerous other considerations that are discussed elsewhere (DBF, [6]) 
but not crucial to our present discussion. 

The smoothing procedure in step 2 furnishes us with a vectorial velocity field. DBF write 
down and use only the radial component of this field. Of course there are reasons to suppose 
the radial component will be more accurately determined than the transverse (^, (f>) components 
since the redshifts are essentially telling us the radial components. However, even if the error 
on the radial component were a factor of ten smaller than on the transverse, it is still possible 
to gain advantage by taking non-radial paths. Figure 1 schematically depicts such a situation. 



3 Errors on the Radial and Transverse Components 



DBF obtain the initial smoothed velocity field, v[f) at arbitrary spatial point r, by minimising 

n 

J2{v.e.r{fi)-UiyW{f,fi) (1) 

8 = 1 

where fi is the position vector of the i^^ galaxy, and there are n galaxies in the catalogue, 
er{fi) is the unit vector in the radial direction at the i^^ galaxy at position fi and W{f,fi) is 
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in accordance with the notation of DBF (appendix A). V is simply the column matrix of 
components of the velocity. Fvidently the inversion of this equation yields all three components 
of the initial smoothed velocity field. 

Since the estimated distances of galaxies are subject to error there will be a corresponding 
error on u(r). Let us write the initial smooth peculiar velocity field obtained from one realisation 
as 

u(r) = u(r) + bv{f) (6) 

where the hat indicates an eshmatec? velocity. DBF's (appendix A) linear error analysis in which 
they derive a bias and variance on the radial component v,^ of u(r) can easily be generalised to 
the vectorial case, although we shall not do this here. The estimated components Or(^) i^4>^^)i 
and ?)e(r) are not statistically independent. If we assume that the bias has been removed from 
?)r(r) so that E[6vr{f)) = 0, we can write the covariance as 



E{6v,{f)6v,{s)) 



Rtj{f,s) 



E denotes expected value. The velocity autocorrelation function, Ri 



(7) 

(r, s), will depend on the 
window function, number density, the dispersion of the distance estimator and also to some 
extent on the input peculiar velocities. Our main purpose in this paper is to demonstrate the 
viability of the method and we shall not attempt to accurately model i?ij(r, s). 



4 Optimal Paths 

The potential, $(r), of the peculiar velocity field is given by the path integral 



v(s).ds 



The error, (5$(r) arises only from the error, 6v[s) in the estimated initial peculiar velocity since 
the path is prescribed. It is therefore given by 

S^{f) = - £ 6iJ{s).ds (9) 

The optimal path for obtaining the potential of the velocity field at position r we shall take to 
be that path for which the variance of the line integral is minimum. 
If we assume that v is unbiased so that E[6v) = , then evidently 

E{6^) = (10) 

The variance of 6^ is given by 

E{S^{rjf = eI f SiJ{s).ds f SiJ{t).dt] = I'' r R,,{s{ji),t{v))x'' {ji)x'' {v)dndv (11) 
\Jo Jo J Jo Jo 

where primes denote derivatives with respect to the arguments. Thus the optimal path will be 
given by 

S r r R^j{s{ii)J{p))x''{ii)x''{p)diidp = (12) 
Jo Jo 

In the case where the autocorrelation function is a delta function, i.e. 

R,,{sJ) = S^{s-t)a,,{s) (13) 

equation (12) simply defines a geodesic on a riemannian space with (Tij as its metric tensor. 
Generally, however, we can expect components of the initial smoothed field at different spatial 
points to be correlated, and the correlation length to be of the same order of magnitude as the 
effective radius of the window function. 

Although equation (12) is very interesting from a mathematical viewpoint, we shall not 
proceed further along those lines. In practice we do not know the form of the autocorrelation 
function, and it could be best approximated by numerical simulations. A more natural way to 
proceed is to use finite element methods, which we now discuss. 



5 Dijkstra's Algorithm 

We wish to calculate the 'best' velocity potential, from which the peculiar velocity field may 
be obtained by taking the gradient. This requires determining the potential $ at regular grid 
points, at least in regions of space where the galaxies are sufficiently dense for the reconstruction 
to be meaningful. Suppose we have N gridpoints at which we wish to evaluate the potential of 
the velocity field. 

Let us assume that the error in going between the a^^ and the b^^ gridpoint is prescribed. 
This 'error length', which we shall call an arcweight, should depend on the number density of 
galaxies in the joint neighborhood of both gridpoints, the distance between the gridpoints and 
the distance of both from the origin. We take the potential to be zero at the first gridpoint (our 
galaxy), and so wish to find the path along which the total error is least. At first sight it might 
appear that this is an NP-complete problem. Luckily this is not the case. Dijkstra's algorithm 
(cf [7]) finds the exact solution to this problem in a number of steps that is bounded above by 
N'^. This algorithm is one of a class of MAXFLOW algorithms used in network theory. 
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Figure 2: Shortest path from to P, indicated by arrows, can be obtained using Dijkstra's 
algorithm. 




Figure 3: Illustration of Dijkstra's algorithm. M is set of nodes for which optimal paths from 
nil have been found. N is the set of all nodes and P = N — M. 



5.1 Description of Dijkstra's Algorithm 



Figure 2 gives a schematic representations of a network of nodes (gridpoints). We write the set 
of nodes as 

N = {ni,n2,n3, ...njv} (14) 

The arcweight between each pair of nodes has been calculated according to some prescription. 
Let the arcweight between node and node nj, be Cab- We shall assume that Cab = c^a, although 
this is not strictly necessary for what follows. A path will be determined by its nodes. Thus 
the path na^na^na^na^-.-naj^ has M nodes and total arcweight or pathlength 

Clearly we may also write 

The problem is to find the path from node rii to every other node that minimises the sum 
of arcweights. We proceed as follows. 

Start from rii. Rename rii as mi (see Figure 3). Suppose at any stage we have a set M of 
nodes to which we have established the minimum pathlength. Let the remaining nodes form 
the set P. Define sm(p) to be the shortest pathlength from mi to p that uses only intermediate 
nodes in M. This must also be the shortest path from mi to p. Now choose from the set of 
nodes P the node p for which the path length is shortest, and add it to the set M and repeat 
until until the set M contains all the nodes. 



6 Application of Dijkstra's algorithm to Potent and Results 

Our main problem in applying this method to POTENT is how to establish the arcweights 
between nodes. The arguments we present now are largely heuristic. Frrors between two nodes 
(gridpoints) we expect to be determined largely by the number density of galaxies in the mutual 
neighbourhood of the two nodes. The higher the galaxy number density the smaller the error. 
Radial components of the initial peculiar velocity field will probably have lower errors than the 
transverse components, but we might expect the two transverse components to have the same 
errors. If we take the gridpoints n and m to be separated by more than the correlation length 
of the autocorrelation function, then we can assume that 6v[m) and 6v[n) are uncorrelated. 
Hence it would be reasonable to take 

E{S^{ f)Y = E{Sv{ab).Ax{ab)Y = ^ E{Sv,{ab)Svj{ab))Ax\ab)Ax^ab) (17) 

where and Ub are consecutive nodes along the path of integration, Ax[ab) the separation be- 
tween these two gridpoints, and 8vi{ab) is the error in the i^^ component of the initial smoothed 
peculiar velocity evaluated at the midpoint of the segment. Since the variance of the distance 
estimator increases with radial distance squared, we shall take the arcweight to also scale with 
. To simplify, we disallow arcs between gridpoints more than three gridlengths away, and 
assume that 

(7re = CTrif, = a^g = and cr^r = k'^crgg = k'^a^^ (18) 

where k is some parameter. 

Thus we shall take Cab to be of the form 



= n"r2((Ar)2 + ^^^^(sin^ e{A^Y + (A^)^) 



(19) 
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Figure 4: Some optimal paths with number density of galaxies indicated by dots. Paths in 
figure (a) have low density weighting and most are almost radial. In figure (b), paths have high 
density weighting and are non-radial. 

where n denotes the number density of galaxies. The value of k in the above equation essentially 
tells one the errors on the transverse components of the initial field compared with the radial. 
This will obviously depend on the window function. It will also depend on the actual peculiar 
velocity field. Very rough simulations indicate that for some simple peculiar velocity field the 
transverse components will be poorly recovered, and consequently k will be large, typically 
between 5 and 10. 

We have taken a spatial distribution of galaxies to be given by [5] combined with [2] data. 
Two peculiar velocity fields are taken corresponding to quiet Hubble fiow and uniform stream- 
ing. Galaxy distances are subjected to distance errors, and the minimum length paths found 
using arcweights of the form (19) to define the line integral (8), and hence the velocity poten- 
tial. Figures 4 and 5 show the optimal paths for two different heuristic arcweight functions, 
and their corresponding rederived velocity fields. 

For both fields it turns out that the optimal paths are almost radial. By insisting on 
arcweights that are heavily dependent on density, and for which A; ~ 1 one can achieve non- 
radial paths. However, in these cases the recovered potential velocity field is noisy and bears 
little resemblance to the input field. 

7 Conclusion and Discussion 

The density inhomogeneities are not great enough to produce optimal paths that are highly 
non-radial when realistic arcweight functions are chosen. Although in regions where the data 
is dense paths deviating from radial do produce potential velocity fields in agreement with 
the radial paths there seems little advantage can be obtained in this way. The reasons for this 
somewhat disappointing result lies primarily with the large radius of the window function which 
tends to smooth out the effect of number density inhomogeneities. For slowly varying peculiar 
velocity fields the transverse components of the initial smoothed peculiar velocity obtained 
from using the window function will be poorly determined, except perhaps near to our own 
galaxy, where in any case the recovery from radial paths is good. Improvements in the derived 
potential velocity could be made by integrating along differentiable curves rather than along 
rectilinear line segments, but probably the gains would not be commensurate with the effort 
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Figure 5: Potential velocities derived from optimal paths obtained for different density weight- 
ings (a) and ratios of transverse to radial errors [k). The solid arrows are the velocities for 
a = and A; = 9. Dashed arrows for a = —4 and k = 5. Galaxy number density is projected 
onto the plane. 



involved. 

One way to improve on the derived potential would be to take the potential at each gridpoint 
averaged over many paths, which could be weighted according to their total arcweight. The 
weakness of POTENT seems to us largely to stem from the arbitrary choice of window function. 
The notion of using an ensemble of paths which are not necessarily radial to obtain an ensemble 
averaged potential is one that can transcend the use of ad hoc window functions, and possibly 
overcome the attendant difficulties of sample gradient and Malmquist bias. We are currently 
investigating this possibility. 
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